//
//  Math.hpp
//  RayTracer
//
//  Created by 俞云康 on 2/15/16.
//  Copyright © 2016 yushroom. All rights reserved.
//

#ifndef Math_h
#define Math_h

#include "RenderFish.hpp"
#include <cmath>
#include <cstdint>
#include <climits>
#include <cfloat>
#include <memory>
#include "Log.hpp"
#include "RNG.hpp"


#ifdef RENDERFISH_PLATFORM_IS_WINDOWS
#define M_E        2.71828182845904523536f
#define M_LOG2E    1.44269504088896340736f
#define M_LOG10E   0.434294481903251827651f
#define M_LN2      0.693147180559945309417f
#define M_LN10     2.30258509299404568402f
#define M_PI       3.14159265358979323846f
#define M_PI_2     1.57079632679489661923f
#define M_PI_4     0.785398163397448309616f
#define M_1_PI     0.318309886183790671538f
#define M_2_PI     0.636619772367581343076f
#define M_2_SQRTPI 1.12837916709551257390f
#define M_SQRT2    1.41421356237309504880f
#define M_SQRT1_2  0.707106781186547524401f
#endif

constexpr float INV_PI = M_1_PI;
constexpr float INV_TWOPI = 1.f / (2 * M_PI);
constexpr float EPS = 1e-3f;
constexpr float EPS2 = EPS * 2;
//constexpr float OneMinusEpsilon = 1 - EPS;

constexpr float MachineEpsilon = std::numeric_limits<float>::epsilon() * 0.5f;

inline bool invalid(float f) {
    return isnan(f) || isinf(f);
}

inline bool is_zero(float f) {
    return fabsf(f) < EPS;
}

inline bool is_one(float f) {
    return fabsf(f-1) < EPS;
}

inline bool isNaN(float f) {
    return isnan(f);
}

inline bool isNaN(int f) {
    return false;
}

inline float radians(float angle) {
    return angle * ( M_PI / 180.0f);
}

inline float angrees(float radian) {
    return radian * ( 180.0f * INV_PI);
}

// Point Declarations
template <typename T>
class Point2 {
public:
    // Point2 Public Methods
    Point2() { x = y = 0; }
    Point2(T xx, T yy) : x(xx), y(yy) { Assert(!HasNaNs()); }
    
    template <typename U>
    explicit Point2(const Point2<U> &p) {
        x = (T)p.x;
        y = (T)p.y;
        Assert(!HasNaNs());
    }
    
//    template <typename U>
//    explicit Point2(const Vector2<U> &p) {
//        x = (T)p.x;
//        y = (T)p.y;
//        Assert(!HasNaNs());
//    }
//    
//    template <typename U>
//    explicit operator Vector2<U>() const {
//        return Vector2<U>(x, y);
//    }
    
#ifndef NDEBUG
    Point2(const Point2<T> &p) {
        Assert(!p.HasNaNs());
        x = p.x;
        y = p.y;
    }
    
    Point2<T> &operator=(const Point2<T> &p) {
        Assert(!p.HasNaNs());
        x = p.x;
        y = p.y;
        return *this;
    }
#endif  // !NDEBUG
    friend std::ostream &operator<<(std::ostream &os, const Point2<T> &p) {
        os << "[" << p.x << ", " << p.y << "]";
        return os;
    }
    
//    Point2<T> operator+(const Vector2<T> &v) const {
//        Assert(!v.HasNaNs());
//        return Point2<T>(x + v.x, y + v.y);
//    }
//    
//    Point2<T> &operator+=(const Vector2<T> &v) {
//        Assert(!v.HasNaNs());
//        x += v.x;
//        y += v.y;
//        return *this;
//    }
//    Vector2<T> operator-(const Point2<T> &p) const {
//        Assert(!p.HasNaNs());
//        return Vector2<T>(x - p.x, y - p.y);
//    }
//    
//    Point2<T> operator-(const Vector2<T> &v) const {
//        Assert(!v.HasNaNs());
//        return Point2<T>(x - v.x, y - v.y);
//    }
//    Point2<T> operator-() const { return Point2<T>(-x, -y); }
//    Point2<T> &operator-=(const Vector2<T> &v) {
//        Assert(!v.HasNaNs());
//        x -= v.x;
//        y -= v.y;
//        return *this;
//    }
    Point2<T> &operator+=(const Point2<T> &p) {
        Assert(!p.HasNaNs());
        x += p.x;
        y += p.y;
        return *this;
    }
    Point2<T> operator+(const Point2<T> &p) const {
        Assert(!p.HasNaNs());
        return Point2<T>(x + p.x, y + p.y);
    }
    template <typename U>
    Point2<T> operator*(U f) const {
        return Point2<T>(f * x, f * y);
    }
    template <typename U>
    Point2<T> &operator*=(U f) {
        x *= f;
        y *= f;
        return *this;
    }
    template <typename U>
    Point2<T> operator/(U f) const {
        float inv = (float)1 / f;
        return Point2<T>(inv * x, inv * y);
    }
    template <typename U>
    Point2<T> &operator/=(U f) {
        float inv = (float)1 / f;
        x *= inv;
        y *= inv;
        return *this;
    }
    T operator[](int i) const {
        Assert(i >= 0 && i <= 1);
        if (i == 0) return x;
        return y;
    }
    
    T &operator[](int i) {
        Assert(i >= 0 && i <= 1);
        if (i == 0) return x;
        return y;
    }
    bool operator==(const Point2<T> &p) const { return x == p.x && y == p.y; }
    bool operator!=(const Point2<T> &p) const { return x != p.x || y != p.y; }
    bool HasNaNs() const { return isNaN(x) || isNaN(y); }
    
    // Point2 Public Data
    T x, y;
};

typedef Point2<int>   Point2i;
typedef Point2<float> Point2f;

struct float2 {
    float x, y;
	float2() : x(0), y(0) {}
	float2(float x, float y) : x(x), y(y) {}
    
    float2 operator+(const float2& v) const {
        return {x+v.x, y+v.y};
    }
    
    float2 operator-(const float2& v) const {
        return {x-v.x, y-v.y};
    }
    float2 operator*(const float f) const {
        return {x*f, y*f};
    }
    
    bool operator<(const float2& v2) const {
        if (x == v2.x) {
            return y < v2.y;
        } else {
            return x < v2.x;
        }
    }
    
    bool operator==(const float2& v2) const {
        return (fabsf(x - v2.x) < 1e-5f && fabsf(y - v2.y) < 1e-5f);
    }

};


struct float3 {
    union {
        struct {
            float x, y, z;
        };
        float m[3];
    };
    
    float3() : x(0), y(0), z(0) { }
    float3(float x, float y, float z) : x(x), y(y), z(z) {
		Assert(valid());
    }
    float3(const float3& b) : x(b.x), y(b.y), z(b.z) { Assert(valid()); }
    
    inline float operator[](int idx) const {
        Assert(idx >= 0 && idx <= 2);
        return m[idx];
    }
    
    inline float& operator[](int idx) {
        Assert(idx >= 0 && idx <= 2);
        return m[idx];
    }

    float3 operator-() const { return {-x, -y, -z}; }
    float3 operator+(const float b)   const { return {x+b, y+b, z+b}; }
    float3 operator-(const float b)   const { return {x-b, y-b, z-b}; }
    float3 operator*(const float b)   const { return {x*b, y*b, z*b}; }
    float3 operator/(const float b)   const { return {x/b, y/b, z/b}; }
    float3 operator+(const float3& b) const { return {x+b.x, y+b.y, z+b.z}; }
    float3 operator-(const float3& b) const { return {x-b.x, y-b.y, z-b.z}; }
    float3 operator*(const float3& b) const { return {x*b.x, y*b.y, z*b.z}; }
    float3 operator/(const float3& b) const { return {x/b.x, y/b.y, z/b.z}; }
    float3& operator+=(const float3& b) { x+=b.x; y+=b.y; z+=b.z; return *this; Assert(valid());}
    float3& operator-=(const float3& b) { x-=b.x; y-=b.y; z-=b.z; return *this; Assert(valid());}
    float3& operator*=(const float3& b) { x*=b.x; y*=b.y; z*=b.z; return *this; Assert(valid());}
    float3& operator+=(const float b) { x+=b; y+=b; z+=b; return *this; Assert(valid());}
    float3& operator-=(const float b) { x-=b; y-=b; z-=b; return *this; Assert(valid());}
    float3& operator*=(const float b) { x*=b; y*=b; z*=b; return *this; Assert(valid());}
    friend float3 operator*(float x, const float3& b) { return float3(x*b.x, x*b.y, x*b.z); }
    //friend float3 operator-(float x, const float3& b) { return {x-b.x, x-b.y, x-b.z}; }
    bool operator<(const float3& v2) const {
        if (x == v2.x) {
            if (y == v2.y) {
                return z < v2.z;
            } else {
                return y < v2.y;
            }
        } else {
            return x < v2.x;
        }
    }
    
    bool operator==(const float3& b) const { return (fabsf(x-b.x)<EPS && fabsf(y-b.y)<EPS && fabsf(z-b.z)<EPS); }
    bool operator!=(const float3& b) const { return (fabsf(x-b.x)>=EPS || fabsf(y-b.y)>=EPS || fabsf(z-b.z)>=EPS); }
    
    bool valid() const { return !(isnan(x)||isnan(y)||isnan(z)||isinf(x)||isinf(y)||isinf(z)); }
    
    float length() const { return sqrtf(x*x + y*y + z*z); }
    float length_square() const { return x*x + y*y + z*z; }
    float normalize_self() { float l = length(); Assert(l > 0 && !isnan(l) && !isinf(l)); float inv_l = 1.f/l; x*=inv_l, y*=inv_l, z*=inv_l; return l;}
    
    float max_element() const { return (x>=y&&x>=z) ? x : (y>z ? y : z); }
};

//typedef float3 Color;
typedef float3 Normal;
typedef float3 Position;
typedef float3 Direction;

struct float4 {
    union {
        struct {
            float x, y, z, w;
        };
        struct {
            float3 xyz; float ww;
        };
        float m[4];
    };
    
    float4() : x(0), y(0), z(0), w(0) { }
    float4(const float x, const float y, const float z, const float w) : xyz(x, y, z), ww(w) {
        Assert(valid());
    }
    float4(const float3& v3, const float w) : xyz(v3), ww(w) {
        Assert(valid());
    }
    float4(const float4& b) : xyz(b.xyz), ww(b.w) {
        Assert(valid());
    }
    float& operator[](int index) {
        Assert(index >= 0 && index <=3);
        return m[index];
    }
    float operator[](int index) const {
        Assert(index >= 0 && index <=3);
        return m[index];
    }
    
    Position as_point() const {
        return xyz / ww;
    }
    
    Direction as_direction() const {
        return xyz;
    }
    
    bool valid() const { return xyz.valid() && !(isnan(z)||isinf(z)); }
};


inline float Dot(const float3& a, const float3& b) { return a.x*b.x + a.y*b.y + a.z*b.z; }


inline float AbsDot(const float3& a, const float3& b) { return fabsf(Dot(a, b)); };


inline float3 Normalize(const float3& v) {
    float inv_l = 1.0f / sqrtf(v.x*v.x + v.y*v.y + v.z*v.z);
    Assert(inv_l > 0 && !invalid(inv_l));
    return float3(v.x*inv_l, v.y*inv_l, v.z*inv_l);
}


inline float3 Cross(const float3& a, const float3& b) {
    return float3(a.y*b.z-a.z*b.y, a.z*b.x-a.x*b.z, a.x*b.y-a.y*b.x);
}


// reflect is normalized if in is normalized
inline float3 Reflect(const float3& in, const float3& normal) {
    return in - (2*Dot(normal, in)) * normal;
}


inline float Distance(const float3& a, const float3& b) {
    return (a-b).length();
}


inline float DistanceSquared(const float3& a, const float3& b) {
    return (a-b).length_square();
}


inline void CoordinateSystem(const Direction &v1, Direction *v2, Direction *v3) {
    if (fabsf(v1.x) > fabsf(v1.y)) {
        float invLen = 1.f / sqrtf(v1.x*v1.x + v1.z*v1.z);
        *v2 = Direction(-v1.z * invLen, 0.f, v1.x * invLen);
    }
    else {
        float invLen = 1.f / sqrtf(v1.y*v1.y + v1.z*v1.z);
        *v2 = Direction(0.f, v1.z * invLen, -v1.y * invLen);
    }
    *v3 = Cross(v1, *v2);
}


inline float clamp(float x, float min_value, float max_value) {
    Assert(min_value < max_value);
    return x<min_value? min_value : x > max_value ? max_value : x;
}


inline float saturate(float x) { return clamp(x, 0, 1); }


inline int   toInt(float x) { return int(powf(saturate(x), 1/2.2f)*255 + .5f); }


inline float Lerp(float t, float v1, float v2) {
    //return (1.f - t) * v1 + t * v2;
    return v1 + (v2 - v1) * t;
}

inline int floor2Int(float val) {
    return (int)floorf(val);
}

inline int ceil2Int(float val) {
    return (int)ceilf(val);
}

inline int round2Int(float val) {
    return floor2Int(val + 0.5f);
}

inline int float2Int(float val) {
    return int(val);
}

inline int Log2Int(uint32_t v) {
#ifdef RENDERFISH_PLATFORM_IS_WINDOWS
	unsigned long lz = 0;
	if (_BitScanReverse(&lz, v))
		return lz;
	return 0;
#else
    return 31 - __builtin_clz(v);
#endif
}

inline constexpr float gamma(int n) {
    return (n * MachineEpsilon) / (1 - n*MachineEpsilon);
}

#endif /* Math_h */
